Hierarchical Time Series

Encuesta de Nacional de Ocupación y Empleo (ENOE)

Motivo de la migración ¹

Los datos son recopilados de la Encuesta Nacional de Ocupación y Empleo (ENOE), siendo, la fuente principal de información sobre el mercado laboral mexicano al ofrecer datos de manera trimestral de la fuerza de trabajo, ocupación, la informalidad laboral, la subocupación y desocupación.

La ENOE inicia su levantamiento en 2005 y se repite en intervalos regulares de tiempo, captando la información a nivel nacional. Este tipo de encuesta permite captar y conocer de mejor manera las características socioeconómicas de la población mexicana de 15 años y más.

El objetivo del trabajo es implementar y comparar las diferentes series de tiempo de manera jerárquica y permitir el desglose de las causas de migración.

Indicadores

Tomando las preguntas del cuestionario de las cuales son más de trescientas con 150 mil casos a nivel nacional en los diferentes periodos de tiempo.

Se toma como referencia la pregunta cs_ad_mot Motivo de la migración en la cual se desglosa en las siguientes categorías.

\(\bullet\) Trabajo
\(\bullet\) Estudio
\(\bullet\) Se casó o unió
\(\bullet\) Se separó o divorció
\(\bullet\) Problemas de salud
\(\bullet\) Reunirse con la familia
\(\bullet\) Inseguridad pública
\(\bullet\) Falleció
\(\bullet\) Otro Motivo
\(\bullet\) No sabe

De las cuales solo se toman 5 casos para el propósito del presente trabajo.

Base de datos

T.Trabajo <- read.xlsx(paste0(here::here(), "/Base de datos/Resultados ENOE.xlsx"),
                       sheet = "Trabajo", 
                       colNames = TRUE, detectDates = TRUE)
T.Estudio <- read.xlsx(paste0(here::here(), "/Base de datos/Resultados ENOE.xlsx"), 
                       sheet = "Estudio",
                       colNames = TRUE, 
                       detectDates = TRUE)
T.Union <- read.xlsx(paste0(here::here(), "/Base de datos/Resultados ENOE.xlsx"),
                     sheet = "Union",
                     colNames = TRUE, 
                     detectDates = TRUE)
T.Divorcio <- read.xlsx(paste0(here::here(), "/Base de datos/Resultados ENOE.xlsx"), 
                        sheet = "Divorcio", 
                        colNames = TRUE, 
                        detectDates = TRUE)
T.Familia <- read.xlsx(paste0(here::here(), "/Base de datos/Resultados ENOE.xlsx"), 
                       sheet = "Reunirse con un familiar", 
                       colNames = TRUE, 
                       detectDates = TRUE)

Se divide entre la población, debido a que si toman los absolutos de las personas que cambiaron de residencia por alguno de los diferentes motivos de ausencia. Estos tienden a ser mayores en algunas ciudades, debido se tiene mayor densidad de población.

Permitiendo así la comparabilidad entre ciudades por una razón de cada 1000 personas.

#Población Total  
T.Poblacion <- read.xlsx(paste0(here::here(), "/Base de datos/Resultados ENOE.xlsx"), 
                         sheet = "Población", 
                         colNames = TRUE, 
                         detectDates = TRUE)
# Vector de tiempo
Periodo <- T.Poblacion$Periodo #Se guarda el vector tiempo

#Creamos una función personalizada para dividir los elementos desde la segunda columna en adelante
divide_columns <- function(df1, df2) {
                    df1_cols <- df1[, 2:ncol(df1)]
                    df2_cols <- df2[, 2:ncol(df2)]
                    
                    # Aplicamos la división usando map2
                    result_cols <- map2_dfc(df1_cols, df2_cols, ~ .x / .y * 1000)
                    
                    # Reconstruimos el data.frame con la primera columna de df1 y las columnas resultantes
                    result <- bind_cols(df1[, 1, drop = FALSE], result_cols)
                    result
}

tablas <- ls(pattern = "T.")
for(i in 1:6){
  assign(paste0(tablas[i]), divide_columns(get(paste0(tablas[i])), T.Poblacion))
}

# Todo en un data.frame
mydata <- do.call(cbind.data.frame, list(T.Trabajo,
                                         T.Estudio %>% select(-c("Periodo")) ,
                                         T.Union %>% select(-c("Periodo")) ,
                                         T.Divorcio %>% select(-c("Periodo")) ,
                                         T.Familia %>% select(-c("Periodo"))))

El tipo de agregación jerárquico es geográfico donde se agrupan las 32 ciudades autorepresentadas en regiones.

Se consideraron 5 regiones del país para que los cálculos no fueran tan extensos a la hora de interpretarlos.

Estructura jerárquica

Cada uno de estos grupos de desglosan en categorías que están anidadas dentro de las categorías de grupos más grandes como 32 ciudades autorepresentadas y 5 regiones, por lo que se hace una recolección de 160 series de tiempo que siguen una estructura de agregación jerárquica.

Esquema Jerárquico

El número de series en el nivel inferior es de 160 series, por lo que las observaciones en el nivel último sumaran a las observaciones del nivel anterior y así sucesivamente.

Los 5 motivos de la ausencia se transforman en series de tiempo, para analizarlos cada uno posteriormente.

ts.Trabajo <- ts(T.Trabajo %>% select(., c(2:length(.))), start = 2005, end = 2019, frequency = 4)
ts.Estudio <- ts(T.Estudio %>% select(., c(2:length(.))), start = 2005, end = 2019, frequency = 4)
ts.Union <- ts(T.Union %>% select(., c(2:length(.))), start = 2005, end = 2019, frequency = 4)
ts.Divorcio <- ts(T.Divorcio %>% select(., c(2:length(.))), start = 2005, end = 2019, frequency = 4)
ts.Familia <- ts(T.Familia %>% select(., c(2:length(.))), start = 2005, end = 2019, frequency = 4)

Motivo de la ausencia: Trabajo

Se analiza cómo caso particular el motivo de la ausencia por trabajo, para empezar familiarizarse con la jerarquización del modelo y posteriormente agregar de manera general los demás motivos expuestos al inicio.

Utilizando la función dygraph permite observar de manera iterativa multiple series de tiempo.

G.Trabajo <- dygraph(ts.Trabajo, main = "Trabajo",  ylab = "Migrantes", xlab = "Periodo") %>%
              dyRangeSelector() %>%
               dyLegend(width = 650) %>%
                dyOptions(colors = RColorBrewer::brewer.pal(32, "Set2")) %>%
                 dyHighlight(highlightSeriesOpts = list(strokeWidth = 2))
G.Trabajo

A continuación, se presenta de manera desagregada un análisis descriptivo de las 32 ciudades autorepresentadas.

Por causas de trabajo la ciudad de Acapulco muestra una media de 4.62 personas que emigraron por trabajo de cada 1000 personas.

Serie de tiempo jerárquica

Se utiliza la función hts para crear una serie de tiempo jerárquica.

Donde se utilizan los datos de último nivel y se utiliza el argumento characters donde los primeros dos caracteres corresponden al primer nivel (Región), los siguientes dos corresponden al segundo nivel (Ciudades) y como último los siete caracteres corresponden al motivo de la ausencia (Trabajo).

#Cambiamos los nombres de las columnas 
Regiones <- c(rep("CE", 7),
              rep("NE", 5), 
              rep("NW", 5), 
              rep("WE", 8),
              rep("SO", 7))
Ciudades <- c("11",     "24",   "01",  "32",    "04",   "14",   "29",   
              "21",     "03",   "06",   "15",   "09",   
              "08",     "19",   "20",   "30",   "18",   
              "12",     "27",   "02",   "05",   "13",   "28",   "22",   "26",   
              "31",     "23",   "07",   "25",   "17",   "10",   "16")
Motivo <- c(rep("Trabajo", 32))

nombres <- paste0(Regiones, Ciudades, Motivo) ## Largo de 11 /Regiones=2,Ciudades=2,Motivo=7
colnames(ts.Trabajo) <- nombres

#Nodos 
Modelo1 <- hts(ts.Trabajo, nodes = list(32), characters = c(2, 2, 7))

## Cambiamos las etiquetas 
Modelo1$labels$`Level 1` <- c("Centro", "Noreste", "Noroeste", "Sureste", "Occidente")

class(Modelo1)
summary(Modelo1)
#p<-smatrix(Modelo1) #Resumen de la matriz de hierarquical time series
#q<-allts(Modelo1) #Matrix de todos los niveles

Nivel 1: Región

Modelo1 %>% 
 aggts(level = 1) %>%
  autoplot(size = 1) + 
   theme_classic() + 
    theme(plot.title = element_text(size = 20),
          plot.subtitle = element_text(size = 12),
          legend.text = element_text(size = 8),
          legend.key.size = unit(0.5, "lines")) +
     scale_color_viridis_d() + 
      scale_color_manual(values = RColorBrewer::brewer.pal(5, "Dark2")) + 
       scale_x_continuous(breaks = seq(2005, 2019, by = 3)) +
        labs(title = "Motivo de la ausencia 2005-2019",
             subtitle = "Trabajo",
             y = "Rate",
             x =  "Year",
             color = "Series") 

La gráfica superior muestra la tasa de migración por trabajo total por región. Donde se puede apreciar que la región Noreste y Noroeste presentan una menor moivilidad con respecto a las demás.

Nivel 2: Ciudad autorepresentada

A continuación se muetran los datos desglosados por las 32 ciudades.

Modelo1 %>% 
 aggts(level = 2) %>%
  autoplot(size = 0.5) + 
   theme_classic() + 
    theme(plot.title = element_text(size = 20, family = "Century Gothic"),
          plot.subtitle = element_text(size = 12, family = "Century Gothic"),
          legend.text = element_text(size = 8, family = "Century Gothic"),
          legend.spacing.x = unit(0.1, "cm"),
          legend.key.size = unit(0.5, "lines"),
          legend.position = "bottom") +
     scale_color_manual(values = colorRampPalette(brewer.pal(8, "Dark2"))(60)) + 
      guides(col = guide_legend(ncol = 15)) +
       labs(title = "Motivo de la ausencia 2005-2019",
            subtitle = "Trabajo",
            y = "Rate",
            x = "Year",
            color = "Series") 

Si bien la gráfica anterior no permite su interpretabilidad debido a que todo el ensamble de las series están muy juntas.

Haciendo uso de la función aggts extrae las series temporales de un objeto hts para cualquier nivel de desagregación.

Para este caso, se están trabajando con dos niveles de desagregación (Regiones y Ciudades) y bien un Nivel Cero (Total).

groups <- aggts(Modelo1, level = 2)

La siguiente gráfica muestra las series de tiempo del nivel inferior, es decir, las tasas de migración por trabajo para cada una de las ciudades en sus respectivas regiones.

Ayudando así a visualizar de manera individual las series.

tibble::as_tibble(groups) %>%
 tidyr::gather(Series) %>%
  mutate(Date = rep(time(groups), NCOL(groups)),
         Group = stringr::str_extract(Series, "([A-Za-z ]*)")) %>%
   ggplot(aes(x = Date, y = value, group = Series, colour = Series)) +
    geom_line() + 
     theme_classic() +
      theme(plot.title = element_text(size = 20, family = "Century Gothic"),
            plot.subtitle = element_text(size = 12, family = "Century Gothic"),
            axis.text.x = element_text(angle = 90, hjust = 1, family = "Century Gothic"),
            legend.text = element_text(size = 8, family = "Century Gothic"),
            legend.spacing.x = unit(0.1, "cm"),
            legend.key.size = unit(0.5, "lines"),
            legend.position = "bottom") +
       scale_color_viridis_d() + 
       #scale_color_manual(values = colorRampPalette(brewer.pal(8, "Set2"))(33)) + 
        scale_x_continuous(breaks = seq(1980, 2015, by = 5)) +
         scale_y_continuous(labels = scales::comma) +
          guides(col = guide_legend(ncol = 15)) +
           labs(title = "Motivo de la ausencia 2005-2019",
                subtitle = "Trabajo",
                y = "Rate",
                x = "Year",
                color = "Grupos") +
            facet_wrap(. ~ Group)

Forecasting: Motivo de trabajo

Utilizando la función forecast() del paquete hts.

Se presentan tres opciones integradas para producir pronósticos usando el parámetro fmethod:

\(\bullet\) ETS Exponential Smoothing,
\(\bullet\) Modelos ARIMA \(\bullet\) Caminatas aleatorias.

Enfoques ²

\(\bullet\) Enfoque de abajo hacia arriba “bottom-up”(method= "bu"): Pronostica del nivel más bajo de la jerarquía, es decir, los motivos de ausencia y luego va agregando los resultados a la jerárquía generar el pronostico del nivel superior al último.

\(\bullet\) Enfoque de arriba hacia abajo “top-down” (method= "tdfp"): Pronostica en la jerarquía del nivel más alto, es decir, por regiones y luego va desglosando los resultados en la jerarquía.

\(\bullet\) Enfoque intermedio **“middle-out” (method= "mo"): Combina enfoques ascendentes y descendentes. Primero, se elige un “nivel medio” y se generan pronósticos para todas las series en este nivel. Para las series por encima del nivel medio, se generan pronósticos coherentes utilizando el enfoque de abajo hacia arriba agregando los pronósticos de “nivel medio” hacia arriba. Para las series por debajo del “nivel medio”, se generan pronósticos coherentes utilizando un enfoque de arriba hacia abajo al desglosar los pronósticos de “nivel medio” hacia abajo.

#h<<-forecast horizon
#method= "mo"<<-"middle-out"
f.modelo1 <- forecast(Modelo1, h = 4, method = "tdfp", fmethod = "ets", parallel = TRUE, keep.fitted = TRUE)

class(f.modelo1)
#summary(f.modelo1)

Nivel 0 al 2

fcst1 <- aggts(f.modelo1, levels = 0:2)
groups <- aggts(Modelo1, levels = 0:2)

autoplot(fcst1, size = 0.5) +
 autolayer(groups) +
  geom_vline(xintercept = 2019, color = "#A8ABD7", linetype= "dashed") +
   theme_classic() +
    theme(plot.title = element_text(size = 18, family = "Century Gothic"),
          plot.subtitle = element_text(size = 14, family = "Century Gothic"),
          legend.text = element_text(size = 7, family = "Century Gothic"),
          legend.key.size = unit(0.5, "lines"),
          legend.position = "bottom") +
     scale_color_viridis_d() + 
      scale_x_continuous(breaks = seq(2005, 2025, by = 2)) +
       scale_y_continuous(labels = scales::comma) +
        guides(col = guide_legend(ncol = 15))+
         labs(title = "Motivo de la ausencia 2005-2019",
              subtitle = "Trabajo",
              y = "Rate",
              x = "Year",
              color = "Series") 

tabla <- ts(rbind(groups, fcst1), start = start(groups), frequency = 4)

Nivel Total: Por motivo de trabajo

#http://www.sthda.com/english/wiki/ggplot2-line-types-how-to-change-line-types-of-a-graph-in-r-software
autoplot(tabla[, "Total"],colour =  "#1720B7", size = 1.2, alpha = 0.6, linetype = "dashed") +
 geom_vline(xintercept = 2019.5, color = "#A8ABD7", linetype = "dashed") +
  theme_classic() +
   scale_x_continuous(breaks = seq(2005, 2025, by = 2)) +
    labs(title = "Motivo de la ausencia 2005-2019",
         subtitle = "Trabajo",
         y = "Rate",
         x =  "Year",
         color = "Series") 

as_tibble(tabla[,-1]) %>%
 tidyr::gather(Series) %>%
  mutate(Date = rep(time(tabla), NCOL(tabla)-1),
         Group = str_extract(Series, "([A-Za-z ]*)")) %>%
   ggplot(aes(x = Date, y = value, group = Series, colour = Series)) +
    geom_line() +
     geom_vline(xintercept = 2017, color = "#A8ABD7", linetype = "dashed") +
      theme_classic() + 
       theme(plot.title = element_text(size = 18, family = "Century Gothic"),
             plot.subtitle = element_text(size = 14, family = "Century Gothic"),
             axis.text.x = element_text(angle = 90, hjust = 1, family = "Century Gothic"),
             legend.position = "bottom",
             legend.text = element_text(size = 8, family = "Century Gothic"),
             legend.key.size = unit(0.5, "lines")) +
        scale_color_viridis_d() +
         scale_x_continuous(breaks = seq(1985, 2025, by = 10)) +
          guides(col = guide_legend(ncol = 15))+
           labs(title = "Motivo de la ausencia 2005-2019",
                subtitle = "Trabajo",
                y = "Rate",
                x = "Year",
                color = "Series") +
            facet_wrap(.~Group) 

f.modelo1 %>% 
 aggts(levels = 0:2) %>%
  autoplot(facet = FALSE) + 
   theme_classic() +
    theme(plot.title = element_text(size = 18, family = "Century Gothic"),
          plot.subtitle = element_text(size = 14, family = "Century Gothic"),
          legend.text = element_text(size = 7, family = "Century Gothic"),
          legend.key.size = unit(0.5, "lines")) +
     scale_color_viridis_d() +
      labs(title = "Motivo de la ausencia 2005-2019",
           subtitle = "Trabajo",
           y = "Rate",
           x = "Year",
           color = "Series")

Motivos de ausencia en general

Utilizando los 5 motivos de ausencia

Serie de tiempo

ts.mydata <- ts(mydata %>% select(., c(2:ncol(.))), start = 2005,end = 2019, frequency = 4)

G.Trabajo<- dygraph(ts.Trabajo, main = "Trabajo", ylab = "Migrantes", xlab = "Periodo") %>%
             dyRangeSelector() %>%
              dyLegend(width = 650) %>%
               dyOptions(colors = RColorBrewer::brewer.pal(32, "Set2")) %>%
                dyHighlight(highlightSeriesOpts = list(strokeWidth = 2))
G.Trabajo

G.Estudio <- dygraph(ts.Estudio, main = "Estudio", ylab = "Migrantes", xlab = "Periodo") %>%
              dyRangeSelector() %>%
               dyLegend(width = 650) %>%
                dyOptions(colors = RColorBrewer::brewer.pal(32, "Set2")) %>%
                 dyHighlight(highlightSeriesOpts = list(strokeWidth = 2))
G.Estudio

G.Union <- dygraph(ts.Union, main = "Se unió o casó", ylab = "Migrantes", xlab = "Periodo") %>%
            dyRangeSelector() %>%
             dyLegend(width = 650) %>%
              dyOptions(colors = RColorBrewer::brewer.pal(32, "Set2")) %>%
               dyHighlight(highlightSeriesOpts = list(strokeWidth = 2))
G.Union

G.Divorcio <- dygraph(ts.Divorcio, main = "Se dicorció o separó", ylab = "Migrantes", xlab = "Periodo") %>%
               dyRangeSelector() %>%
                dyLegend(width = 650) %>%
                 dyOptions(colors = RColorBrewer::brewer.pal(32, "Set2")) %>%
                  dyHighlight(highlightSeriesOpts = list(strokeWidth = 2))
G.Divorcio

G.Familiar <- dygraph(ts.Familia, main = "Reunirse con un familiar", ylab = "Migrantes", xlab = "Periodo") %>%
               dyRangeSelector() %>%
                dyLegend(width = 650) %>%
                 dyOptions(colors = RColorBrewer::brewer.pal(32, "Set2")) %>%
                  dyHighlight(highlightSeriesOpts = list(strokeWidth = 2))
G.Familiar

Modelo jerárquico

Se utiliza la función hts() para crear una serie de tiempo jerárquica.

Donde se utilizan los datos de último nivel y se utiliza el argumento characters donde los primeros dos caracteres corresponden al primer nivel (Región), los siguientes dos corresponden al segundo nivel (Ciudades) y como último los tres caracteres corresponden al motivo de la ausencia.

\(\bullet\) TRA= “Trabajo”
\(\bullet\) EST= “Estudio”
\(\bullet\) UNI= “Se casó o unió”
\(\bullet\) DIV= “Se divorció o separó”
\(\bullet\) FAM= “Reunirse con un familiar”

Regiones <- rep(c(rep("CE", 7),
                  rep("NE", 5),
                  rep("NW", 5),
                  rep("WE", 8),
                  rep("SO", 7)), 5)
Ciudades <- rep(c("11",     "24",   "01",   "32",   "04",   "14",   "29",
                  "21",     "03",   "06",   "15",   "09",
                  "08",     "19",   "20",   "30",   "18",   
                  "12",     "27",   "02",   "05",   "13",   "28",   "22",   "26",   
                  "31",     "23",   "07",   "25",   "17",   "10",   "16"), 5)

Motivo <- c(rep("TRA", 32),
            rep("EST", 32),
            rep("UNI", 32),
            rep("DIV", 32),
            rep("FAM", 32))

nombres <- paste0(Regiones, Ciudades, Motivo) #Largo de 7 | Regiones=2|Ciudades=2|Motivo=3
colnames(ts.mydata) <- nombres

nodes <- list(160, c(5, 32, 5)) #160 variables |Regiones=5|Ciudades=32|Motivos=5|  
Modelo2 <- hts(ts.mydata, nodes = nodes,characters = c(2, 2, 3))

#Cambiamos los labels 
Modelo2$labels$`Level 1` <- c("Centro", "Noreste", "Noroeste", "Sureste", "Occidente")

Nivel 1: Región

Modelo2 %>% 
 aggts(level = 1) %>%
  autoplot(size = 1) + 
   theme_classic() + 
    theme(plot.title = element_text(size = 20),
          plot.subtitle = element_text(size = 12),
          legend.text = element_text(size = 8),
          legend.key.size = unit(0.2, "cm")) +
     scale_color_viridis_d() + 
      scale_color_manual(values = RColorBrewer::brewer.pal(5, "Dark2")) + 
       scale_x_continuous(breaks = seq(2005, 2019, by = 2)) +
        labs(title = "Motivo de la ausencia 2005-2019",
             subtitle = "Trabajo / Estudio / Se casó o unió / Divorció o separó / Reunirse con un familiar",
             y = "Rate",
             x = "Year",
             color = "Series") 

Nivel 2: Cuidad autorepresentada

Modelo2 %>% 
 aggts(level = 2) %>%
  autoplot(size = 0.5) + 
   theme_classic() + 
    theme(plot.title = element_text(size = 20),
          plot.subtitle = element_text(size = 12),
          legend.text = element_text(size = 8),
          legend.key.size = unit(0.5, "lines"),
          legend.position = "bottom") +
     scale_color_manual(values = colorRampPalette(brewer.pal(8, "Dark2"))(60)) + 
      guides(col = guide_legend(ncol = 15))+
       labs(title = "Motivo de la ausencia 2005-2019",
            subtitle = "Trabajo / Estudio / Se casó o unió / Divorció o separó / Reunirse con un familiar",
            y = "Rate",
            x = "Year",
            color = "Series") 

Si bien la gráfica anterior no permite su interpretabilidad debido a que todo el ensamble de las series están muy juntas.

Haciendo uso de la función aggts extrae las series temporales de un objeto hts para cualquier nivel de desagregación.

groups <- aggts(Modelo2, level = 2)

tibble::as_tibble(groups) %>%
 tidyr::gather(Series) %>%
  mutate(Date = rep(time(groups), NCOL(groups)),
         Group =stringr::str_extract(Series, "([A-Za-z ]*)")) %>%
   ggplot(aes(x = Date, y = value, group = Series, colour = Series)) +
    geom_line() + 
     theme_classic() +
      theme(plot.title = element_text(size = 20),
            plot.subtitle = element_text(size = 12),
            axis.text.x = element_text(angle = 90, hjust = 1),
            legend.text = element_text(size = 8),
            legend.key.size = unit(0.5, "lines"),
            legend.position = "bottom") +
       scale_color_viridis_d() +
       #scale_color_manual(values = colorRampPalette(brewer.pal(8, "Dark2"))(33)) + 
        scale_x_continuous(breaks = seq(1980, 2015, by = 5)) +
         scale_y_continuous(labels = scales::comma) +
          guides(col = guide_legend(ncol = 15))+
           labs(title = "Motivo de la ausencia 2005-2019",
                subtitle = "Trabajo / Estudio / Se casó o unió / Divorció o separó / Reunirse con un familiar",
                y = "Rate",
                x = "Year",
                color = "Series") +
            facet_wrap(. ~ Group) 

Nivel 3: Motivo de la ausencia

En modelo general, se están trabajando con tres niveles de desagregación (Motivo de trabajo, Regiones y Ciudades) y bien un Nivel Cero (Total).

groups <- aggts(Modelo2, level = 3)

#https://stringr.tidyverse.org/reference/str_locate.html
tibble::as_tibble(groups) %>%
 tidyr::gather(Series) %>% 
  mutate(Date = rep(time(groups), NCOL(groups)),
         Group =stringr::str_extract(Series, "([A-Za-z ]*)"),
         Motivo=stringr::str_sub(Series,5, 7)) %>%
   ggplot(aes(x = Date, y = value, group = Series, colour = Series)) +
    geom_line() + 
     theme_classic() +
      theme(plot.title = element_text(size = 20),
            plot.subtitle = element_text(size = 12),
            axis.text.x = element_text(angle = 90, hjust = 1),
            legend.text = element_text(size = 8),
            legend.key.width=unit(0.2, "cm"),
            legend.key.height = unit(0, "cm"),
            legend.spacing.x = unit(0, "cm"),
            legend.key.size = unit(0.5, "lines"),
            legend.position = "bottom") +
       scale_color_viridis_d() +
       #scale_color_manual(values = colorRampPalette(brewer.pal(8, "Dark2"))(160)) + 
        scale_x_continuous(breaks = seq(1980, 2015, by = 5)) +
         scale_y_continuous(labels = scales::comma) +
          guides(col = guide_legend(ncol = 15))+
           labs(title = "Motivo de la ausencia 2005-2019",
                subtitle = "Trabajo / Estudio / Se casó o unió / Divorció o separó / Reunirse con un familiar",
                y = "Rate",
                x = "Year",
                color = "Series") +
            facet_wrap(. ~ Group + Motivo) 

Analizando la gráfica anterior el motivo de ausencia para “Reunirse con un familiar” es la que presenta un mayor número de casos migratorios” con respecto a los otros cuatro casos.

#print(Modelo2)
#smatrix(Modelo2)
#allts(Modelo2)

Forecasting: Motivo de la ausencia en general

Se espera que los pronósticos sean consistentes con la estructura de agregación de las series de tiempo al agruparlas.

#h=forescast horizon
#method<-"mo", "bu", "tdfp"
f.modelo2 <- forecast(Modelo2, h = 10, method = "mo", level = 2, fmethod = "arima", parallel = TRUE, keep.fitted = TRUE)

summary(f.modelo2)

Nivel 0 al 2: Total / Región / Ciudad

fcst2 <- aggts(f.modelo2, levels = 0:2)
groups <- aggts(Modelo2, levels = 0:2)

autoplot(fcst2, size = 0.5) +
 autolayer(groups) +
  geom_vline(xintercept = 2019,color = "#A8ABD7", linetype = "dashed") +
   theme_classic() +
    theme(plot.title = element_text(size = 20),
          plot.subtitle = element_text(size = 12),
          legend.text = element_text(size = 7),
          legend.key.size = unit(0.5, "lines"),
          legend.position = "bottom") +
     scale_color_viridis_d() + 
      scale_x_continuous(breaks = seq(2005, 2025, by = 2)) +
       scale_y_continuous(labels = scales::comma) +
        guides(col = guide_legend(ncol = 15))+
         labs(title = "Motivo de la ausencia 2005-2019",
              subtitle = "Trabajo / Estudio / Se casó o unió / Divorció o separó / Reunirse con un familiar",
              y = "Rate",
              x = "Year",
              color = "Series") 

tabla <- ts(rbind(groups, fcst2),
                  start = start(groups), frequency = 4) #Frecuencia al año

Nivel 0: Total de casos de migración

#http://www.sthda.com/english/wiki/ggplot2-line-types-how-to-change-line-types-of-a-graph-in-r-software
autoplot(tabla[, "Total"], colour = "#1720B7", size = 1.2, alpha = 0.6, linetype = "dashed") +      
 geom_vline(xintercept=2019.5,color = "#A8ABD7",linetype= "dashed") +
  theme_classic() +
   theme(plot.title = element_text(size = 20),
         plot.subtitle = element_text(size = 12)) + 
     scale_x_continuous(breaks = seq(2005, 2025, by = 2)) +
     labs(title = "Motivo de la ausencia 2005-2019",
          subtitle = "Trabajo / Estudio / Se casó o unió / Divorció o separó / Reunirse con un familiar",
          y = "Rate",
          x = "Year",
          color = "Series") 

Nivel 1 a 2: Nivel región y nivel ciudad

as_tibble(tabla[,-1]) %>%
 tidyr::gather(Series) %>%
  mutate(Date = rep(time(tabla), NCOL(tabla)-1),
         Group = str_extract(Series, "([A-Za-z ]*)")) %>%
   ggplot(aes(x = Date, y = value, group = Series, colour = Series)) +
    geom_line() +
     geom_vline(xintercept = 2019,color = "#A8ABD7", linetype = "dashed") +
      theme_classic() + 
       theme(plot.title = element_text(size = 20),
             plot.subtitle = element_text(size = 12),
             axis.text.x = element_text(angle = 90, hjust = 1),
             legend.text = element_text(size = 8),
             legend.key.width = unit(0.2, "cm"),
             legend.key.height = unit(0, "cm"),
             legend.spacing.x = unit(0.1, "cm"),
             legend.key.size = unit(0.5, "lines"),
             legend.position = "bottom") +
        scale_color_viridis_d() +
         scale_x_continuous(breaks = seq(2005, 2025,by = 5)) +
          guides(col = guide_legend(ncol = 15)) +
           labs(title = "Motivo de la ausencia 2005-2019",
                subtitle = "Trabajo / Estudio / Se casó o unió / Divorció o separó / Reunirse con un familiar",
                y = "Rate",
                x = "Year",
                color = "Series") +
            facet_wrap(. ~ Group) 

Nivel 3: Motivos de ausencia en general

fcst3 <- aggts(f.modelo2, levels = 3)
groups <- aggts(Modelo2, levels = 3)
tabla <- ts(rbind(groups, fcst3), start = start(groups), frequency = 4) #Frecuencia al año

as_tibble(tabla[,-1]) %>%
 tidyr::gather(Series) %>%
  mutate(Date = rep(time(tabla), NCOL(tabla) - 1),
         Group = str_extract(Series, "([A-Za-z ]*)"),
         Motivo = stringr::str_sub(Series, 5, 7)) %>%
   ggplot(aes(x = Date, y = value, group = Series, colour = Series)) +
    geom_line() +
     geom_vline(xintercept = 2019, color = "#A8ABD7", linetype= "dashed") +
      theme_classic() + 
       theme(plot.title = element_text(size = 20),
             plot.subtitle = element_text(size = 12),
             axis.text.x = element_text(angle=90, hjust=1),
             legend.text = element_text(size = 6),
             legend.key.width=unit(0.2, "cm"),
             legend.key.height = unit(0, "cm"),
             legend.spacing.x = unit(0.1, "cm"),
             legend.key.size = unit(0.01, "lines"),
             legend.position = "bottom") +
        scale_color_viridis_d() +
         scale_x_continuous(breaks = seq(2005, 2025,by = 5)) +
          guides(col = guide_legend(ncol = 15)) +
           labs(title = "Motivo de la ausencia 2005-2019",
                subtitle = "Trabajo / Estudio / Se casó o unió / Divorció o separó / Reunirse con un familiar",
                y = "Rate",
                x = "Year",
                color = "Series") +
            facet_wrap(.~ Group + Motivo)

Referencias

Encuesta Nacional de Ocupación y Empleo (ENOE), población de 15 años y más de edad. (n.d.). Retrieved March 30, 2020, from https://www.inegi.org.mx/programas/enoe/15ymas/

Forecasting Hierarchical Time Series using R - Brillio Data Science - Medium. (n.d.). Retrieved March 30, 2020, from https://medium.com/brillio-data-science/forecasting-hierarchical-time-series-using-r-598828dba435

Librerías

Librerías que se usaron en el trabajo

Nos ha servido a construir modelos de series de tiempo con datos estructurales para fines de pronósticos.

This work by [Diana Villasana Ocampo]{xmlns:cc= “http://creativecommons.org/ns#” property = “cc:attributionName”} is licensed under a Creative Commons Attribution 4.0 International License.

---

1. (https://www.inegi.org.mx/programas/enoe/15ymas/)↩︎

2. (Forecasting Hierarchical Time Series using R - Brillio Data Science - Medium. (n.d.). Retrieved March 31, 2020, from https://medium.com/brillio-data-science/forecasting-hierarchical-time-series-using-r-598828dba435)↩︎